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The detection of deep reflected S waves on Mars inferred a core size of 1,830 + 40 km 
(ref. 1), requiring light-element contents that are incompatible with experimental 
petrological constraints. This estimate assumes a compositionally homogeneous 


Martian mantle, at odds with recent measurements of anomalously slow propagating 
P waves diffracted along the core-mantle boundary’. An alternative hypothesis is that 
Mars’s mantle is heterogeneous as a consequence of an early magma ocean that 
solidified to form a basal layer enriched in iron and heat-producing elements. Such 
enrichment results in the formation of a molten silicate layer above the core, overlain 
by a partially molten layer’. Here we show that this structure is compatible with all 
geophysical data, notably (1) deep reflected and diffracted mantle seismic phases, 
(2) weak shear attenuation at seismic frequency and (3) Mars’s dissipative nature at 
Phobos tides. The core size in this scenario is 1,650 + 20 km, implying a density of 

6.5 g cm”, 5-8% larger than previous seismic estimates, and can be explained by 
fewer, and less abundant, alloying light elements than previously required, inamounts 
compatible with experimental and cosmochemical constraints. Finally, the layered 
mantle structure requires external sources to generate the magnetic signatures 
recorded in Mars’s crust. 


Mars, like other differentiated terrestrial planets, is composed of an 
iron alloy core overlain by a silicate mantle and a crust. This metal- 
silicate dichotomy probably originates from an early global magma 
ocean stage*» during which heavy iron gravitationally separates from 
the lighter silicates to form a core®. Observations collected by space 
missions have improved our knowledge on Mars’s present-day internal 
structure. Among these, the NASA InSight mission’ that deployed the 
first seismometer on the surface of Mars has been instrumental’. Quake 
recordings’? have determined the layering of the crust", its thickness”, 
mantle structure», the core size of Mars and its composition”, 
The seismic detection of Mars’s core confirmed the geodetic measure- 
ments”! with an inferred core radius of R, = 1,830 + 40 km (ref. 1), 
recently revised to slightly smaller values (R, = 1,780 - 1,810 km)”, and 
implies a relatively low core density of approximately 6.0-6.2 g cm”. 
This core size was determined by the detection of S waves reflected at 
asolid-liquid interface located at the bottom of the solid mantle, and 
ascribed to be the core-mantle boundary (CMB). Like all InSight pre- 
and post-mission structure models of Mars, the mantle was assumed 
to be compositionally homogeneous! "6, 

The presence of alloying light elements can explain alow core density, 
though cosmochemical and experimental constraints do limit this 
explanation. Sulfur (S) alone could produce the required densities, 
but with concentrations significantly above the maximum allowable 
by cosmochemistry”. Hence, additional elements have been invoked to 


reduce the amount of S, with oxygen concentrations around 5 wt% and 
carbon and hydrogen concentrations of approximately 1 wt% (ref. 1). 
However, such concentrations of hydrogen and carbon are significantly 
higher than allowed by experimental observation. Indeed, experimen- 
tal constraints” limit the hydrogen content in the core to 0.15 wt%. In 
addition, sulfur decreases the solubility of carbon in Fe-rich alloys from 
more than 4 wt% (for an S-free core) to below 1 wt% (ref. 29) for a core 
containing 16 wt% of sulfur. Finally, oxygen concentration in molten 
Fe-S alloys is controlled by the sulfur content”, and cannot be larger 
than 4 wt%, a value well below that proposed in ref. 1. 

Since the first seismic determination of Mars’s present-day struc- 
turet? additional events have been recorded by InSight. In par- 
ticular, on Sol (a Martian day, whose duration is 24 h 40 min) 1,000 
after InSight’s landing, aseismic event triggered by a meteorite impact 
(hereafter named S1000a) was detected. The impact was precisely 
located by orbital imaging at 125.9° away from InSight’, and triggered 
seismic signals interpreted as the first-ever observation of P waves dif- 
fracted along Mars’s CMB (hereafter labelled Pdiff)". Among hundreds 
of models published!**+ only a handful of models can marginally fit 
the measured differential arrival times between Pdiff and PP (that is, 
a P wave reflected once at the surface of the planet)’. A proper fit 
requires a significant velocity reduction in the deep mantle, which 
cannot be explained in a homogeneous mantle commonly assumed 
in Mars’s structure models*”’, 
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Fig. 1|Inversion results. Thermochemical evolution and present-day structure 
of Mars. One model among the best 50 is displayed for each inversion set. 

a-f, Without a BML (homogeneous mantle), with n= 6 x 10% Pa s, E*=300 kJ mol”, 
V*=3.8 cm’ mol”. g-I, With a BML (heterogeneous mantle), with nọ=5 x 10° Pas, 
E*=110 kJ mol”, V*= 4.4 cm? mol”. a,g, Evolution of crustal and lithospheric 
thicknesses (including the uppermost mantle thermal boundary layer (TBL)). 
b,h, Evolution of uppermost convective mantle (7,,) and core (T.) temperatures. 
c,i, Present-day temperature profiles and mantle melting curves from ref. 48 
accounting for the influence of iron inthe BML’. d,j, Density profiles. e,k, Shear 
and compressional wave speed profiles. f,1, Raypaths for waves reflected 

at (blue) or diffracted along (red) deep mantle interfaces. Additional raypaths 
for other phases considered for the inversion are shown in grey. m,n, Close-up 
views of the region delineated in the vicinity of the BML in l, showing the P- and 


Additionally, the detection of distant small magnitude events by 
InSight suggests that seismic attenuation is weak, with effective shear 
quality factors in excess of 1,000 (ref. 9). This contrasts with the attenu- 
ation measured at larger periods such as during Phobos’s main tides 
(5h 55 min), with a corresponding global quality factor Q, = 95 +10 
(refs. 34,35) indicative ofa relatively attenuating mantle for tides. This 
behaviour may be explained by a stratified mantle, akin to a soft, deep 
mantle layer, as proposed for the Moon”**. 


S-wave velocity structure (m) and raypath (n) of the P-diffracted wave reflected 
at the CMB (PbdiffPcP). In the homogeneous mantle, S-wave reflection occurs 
at the CMB, while in the heterogeneous mantle, it occurs above the CMB where 
velocity decreases abruptly due to the transition froma partially molten to 
afully/essentially molten state in the BML (dotted curves). Inthe heterogeneous 
mantle, the P-diffracted phase (PbdiffPcP) travels in a molten silicate mantle 
with slower wave speeds compared with those ina solid mantle, significantly 
delaying its travel time. The PbdiffPcP phase results from multiple rays 
diffracted at the top of the fully molten BML before and/or after core reflection, 
which contribute to this seismic phase. The path displayed corresponds to one 
of these contributions, which is the reason why it is not symmetric. However, 
because the seismic model is spherically symmetric, the sum of the contributions 
will result in asymmetric path (Supplementary Fig. 4). 


The validity of the common assumption of ahomogeneous Martian 
mantle is brought into question by seismic, geodetic, cosmochemical 
and experimental observations. In fact, the solidification of Mars’s early 
magma ocean could produce a heterogeneous mantle, leading toa 
deep-rooted silicate layer strongly enriched in iron and heat-producing 
elements (HPEs) just above the core” “, which is also suggested by 
isotopic anomalies measured in Martian meteorites” “*. Such a basal 
mantle layer (hereafter BML) can strongly affect Mars’s thermochemical 
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evolution, and may considerably influence the interpretation of avail- 
able geophysical data*™. Its presence leads to the development of a 
soft, partially molten mantle layer that could explain Mars’s dissipating 
behaviour for tides®. It also implies that the observed S-wave reflections 
occur near the top of the BML and not at the CMB, suggesting a smaller 
and denser core than the recently inferred values, which appears com- 
patible with the recently measured nutation of Mars*. These lines of 
evidence suggesting the existence of a previously unrecognized BML 
motivate the reinterpretation of available data used to constrain the 
interior structure of Mars. To determine planet structures compatible 
with observations, we performed a probabilistic inversion of seismic 
data (Supplementary Information Sections 1 and 3). We considered a 
non-BML inversion set, with a compositionally homogeneous mantle, as 
commonly assumed for Mars, anda second set including a BML above 
the core’. Our inversion is parameterized in terms of quantities that 
influence the thermochemical evolution of the planet composed ofa 
liquid iron core, a silicate mantle (with or without a BML), and an evolv- 
ing lithosphere and crust?”°**. Mantle viscosity controls the thermal 
evolution of the planet and depends ontemperature, 7, and pressure, 
P(for example, ref. 47): 


E* + P(r) “), (1) 


n(r) = Aref ner RT(r) 
where the effective activation energy E* expresses the temperature 
sensitivity, the effective activation volume V* controls the pressure 
dependence, R is the gas constant, no is the reference viscosity and 
Aeris a prefactor. 

Figure 1 underlines the differences between standard (Fig. la—f) and 
BML models (Fig. 1g-1). For ahomogeneous mantle, core cooling is 
monotone, while the uppermost mantle temperature (Tmn) firstincreases 
dueto HPEs, before decreasing again. The hot (approximately 1,760 K) 
initial mantle temperatures favour a rapid crustal growth that depletes 
the mantle in heat sources, leading to mantle cooling (Fig. 1a,b). This 
yields a T,, approximately 100 K hotter than initially (Fig. 1b). The 
present-day temperature profile (Fig. 1c) results in mineralogies and 
associated density (Fig. 1d) and seismic velocity profiles (Fig. le). The 
latter lead to S-wave reflections occurring at the CMB, as in ref. 1, and 
a Pdiff phase propagating along the CMB (Fig. 1f). 

The BML acts as a buffer that reduces heat transfer between the 
core andthe mantle, but also as a heat source’, resulting in substantial 
core heating from approximately 2,160 K to approximately 2,840 K 
(Fig. 1h). Since the BML segregates most of the HPEs, heat sources are 
depleted in the overlying mantle, leading to a rapid mantle cooling 
from approximately 1,820 K to approximately 1,530 K. Yet, the mantleis 
initially about 60 K hotter than in the homogeneous model (Fig. 1g-1). 
This favours shallow mantle melting and low mantle viscosities (equa- 
tion (1)), leading to more efficient cooling, a crustal formation com- 
pleted more than 500 Myr earlier, and a lithospheric growth about two 
times slower (Fig. 1g) than the homogeneous case (Fig. 1a). The BML HPE 
enrichment leads to temperatures exceeding melting curves*’, result- 
ing ina present-day BML fully and partially molten, overlain by a cold 
(7,, approximately 1,530 K) depleted mantle (Fig. 1i). This present-day 
temperature profile leads to a distinct density structure and faster 
velocities in the mantle above the BML but strongly reduced P- and 
S-wave velocities (Vp and V,) inthe hot and molten BML (Fig. 1j-k). The 
zero V,; in the lower part of the BML yields an S-wave reflection occur- 
ring 150 km above the CMB, at the interface where the melt fraction 
becomes large enough to behave as a liquid material*’™ (Fig. 11). In 
this case, there exists a diffracted P wave in the deep mantle, whose 
path differs significantly from the one associated with the non-BML 
case (Fig. 1f). This P wave diffracts along the bottom of the depleted 
mantle just above the BML, then travels down through the BML before 
reflecting back at the CMB towards the surface (Fig. 1I-n). We call 
this core-bouncing diffracted phase PbdiffPcP (whose occurrence is 


714 | Nature | Vol 622 | 26 October 2023 


a b 
3,000 F 3,000 + 
Ẹ Ẹ 
2 Erd 
w 2,500 F m 2500F 
=] =| 
re] xe} 
oO oO 
fia fam 
2,000 F Without a BML 2,000 | With a BML 
m fe ee eee 1 
| l |] l 
| | | ! 
J Hn) ---F____ {l 
1,500 === ee n] 1,500 : 
0 2,500 5,000 7,500 10,000 0 2,500 5,000 7,500 10,000 
Vp or Vs (m s~’) Vp or V; (m s~’) 
c 1,900 - d 
CMB range 
Ẹ 1,800 F Ẹ 
= = —— 
3 s Fully molt 
2 2 ully molten 
3 1700F 8 1,700 F BML range 
ICMB range | i 
1,600 + 1,600 + ad 
J L 1 J 
0 2,500 5,000 7,500 10,000 0 2,500 5,000 7,500 10,000 
Vp or Vg (m s~’) Vp or Vg (m s~’) 
e f 
Œ R, (without a BML) ee 
_ 15+ æ R, (with a BML) — 15 - m Without a BML 
g m RY (with a BML) g m With a BML 
> > 
g 107 | S aot | 
G i S i 
=) Pi 
fom Å [og { | 
D fi K mal 
È Sr J i Ù Sr [i 
fi | i. 1 iz i 
1,600 1,700 1,80 1,900 6,000 6,250 6,500 6,750 


Real (R,) or apparent (R) Average core density, p, (kg m™) 


core radius (km) 
Fig. 2 | Inversion results of Mars’s seismic data with or without a BML. 
a, Seismic profiles for P-wave (red) and S-wave (blue) velocities for 50 models 
chosen among the best 1,000 models without a BML. b, Same as a but fora 
mantle that contains a BML. c,d, Zoom-in of the areas delineated by the dashed 
rectangles in a and b, respectively. The maximum or minimum depth ranges of 
three distinct seismic regions and interfaces are marked by vertical arrows: the 
CMB range, the range for the interface between the mushy layer and the fully 
molten BML, and the fully/essentially molten BML region. e, Histograms for 
real and apparent core radii for the best 1,000 models. The core radius is 
considerably smaller when a BML is present but the apparent core radii (that is, 
the radius of liquid iron alloy plus the thickness of the fully molten silicate layer) 
are similar in both cases. f, Histograms of Mars’s core density for the best 1,000 
models. The smaller core size in the heterogeneous-mantle case leads toa denser 
core compared with the homogeneous-mantle case. Panels a-d contain a smaller 
number of models compared with e and fto allow for a clear visualization of the 
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demonstrated in Supplementary Information Section 2). The PbdiffPcP 
wave travels within the molten silicate layer (Fig. 1m,n), thereby increas- 
ing its travel time significantly (Fig. 1k). 

Further differences emerge between BML and non-BML sets 
(Extended Data Table 1). The present-day mantle including a BML is 
approximately 300 K colder because the mantle is, on average, initially 
70 K hotter and about 30 times less viscous than non-BML models. 
In addition, BML models have £* and V* about three times smaller and 
about 60% larger, respectively, than those associated with non-BML 
models. The combination of a cold mantle with a small £* for the BML 
set leads to uppermost present-day mantle viscosities comparable 
to those in homogeneous models (n = 10” Pa s). The differences in V* 
imply that the convecting mantle (regions displayed in orange in Fig. 1) 
above the BML becomes considerably more viscous. Finally, the initially 
hotter mantle, but now colder in the present-day, results in a thinner 
lithosphere for the BML set. 
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Fig. 3 | Datafit for the differential travel times, considering a mantle 
without and with a BML for the best 1,000 models. The results are displayed 
interms of probability density functions (PDFs). Blue and red colours show 
small and large probabilities, respectively. a,b, Body wave differential times 
with respect to P waves as a function of ts- tp. Markers and error bars represent 
observed differential arrival times and their uncertainty, respectively. 


Core quantities also differ between the two inversion sets. For the 
BML set R, is approximately 170 km smaller (R, = 1,650 + 20 km ver- 
sus R, = 1,820 + 10 km) and therefore denser (average core density 
Pp. = 6,470 + 60 kg m”? versus p. = 6,130 + 30 kg m°) than the non-BML 
set to satisfy constraints on the mass of Mars. The CMB temperature 
is approximately 700 K hotter for the BML set. 

While the seismic structures for models with homogeneous mantles 
simply consist ofa liquid core with V; = O ms ‘overlain by a solid silicate 
envelope (Fig. 2a,d), the seismic structures for the BML models are more 
complicated. The increasing temperature with depth (Fig. 1i) causes an 
increase in melt fraction with depth in the BML. Therefore, the layer is 


C-f, top— tpaire (C), Cop — Copairrpep (d) and tsks~ tpp (€,f) as a function of tss~ tpp for event 
$1000a. The observed differential travel time measurements are displayed 
with black lines. The pink bands indicate uncertainties on the measurements. 
Panels a,c,e correspond to the inversion set without a BML, and b,d,f to the 
inversion set witha BML. 


subdivided into a mushy part (that is, with a melt fraction below 65% 
(ref. 50), implying strongly reduced but non-zero V;) overlying a fully/ 
essentially molten part (that is, with V; = 0 ms‘) in various propor- 
tions (Fig. 2d). The core radius plus the thickness of the fully molten 
silicate layer (that is, the apparent core radius, R, = 1,780 + 20 km) is 
comparable, yet approximately 40 km smaller than the core size for the 
non-BML inversion set. However, R, leans towards recent re-estimation 
of R.=1,780 - 1,810 km (median range) using a larger dataset”. 
Constraints from S waves reflected off a deep solid-liquid inter- 
face yield a core on average 170 km smaller than the homogeneous 
case (Fig. 2e). Yet, the BML apparent radii are comparable to R, in the 
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Fig. 4| Attenuation at seismic and Phobos’s tide periods deduced from the 
100 best models. a,b, Present-day shear attenuation profiles at 1 Hz for the 
inversion set without a BML (a) and witha BML (b) for various frequency 
dependence (a,) of the shear quality factor. c, Apparent shear quality factor for 
the S-reflected wave, Q,.s, at 1 Hz and for event S1222a as a function of the 
global quality factor at Phobos tidal frequency (5 h 55 min). Qses was computed 
by integrating 1/Q, for each model along the ScS raypath, determined using the 
TauP toolkit®°. Models without BML are displayed in green tones and models 
witha BML areas shown in purple tones depending on the value of a,. The 
acceptable ranges for Q,., and Q,, are displayed in pink. d, Q.,,at 1 Hz and for 
event S1222a as a function of Rọ, the ratio of reflected S to direct S amplitudes 
Ages/As between events S0185a and S1222a at 1 Hz that samples different solid 
mantle depths. The colour and symbol coding are identical to that in c. 


homogeneous cases, due to travel time constraints from deep reflected 
S waves. As seen above, the smaller core in the BML set has a higher 
density than homogeneous models (Fig. 2f). In agreement with ref. 45 
this BML density structure is compatible with constraints from Mars’s 
core nutation (Supplementary Information Section 6), and models with 
large R, display the smallest misfit to observations”. 

While both inversion outputs generally exhibit good fit to seismic 
data (Fig. 3a,b), the BML models specifically yield a considerably bet- 
ter fit to waves propagating in the lowermost mantle and in the core. 
For a homogeneous mantle, tpp — tpaite falls outside the lo uncertainty 
(Fig. 3c), confirming that the P wave diffracted along the CMB travels 
too fast relative to the PP-phase’. By contrast, the molten BML layer 
significantly delays the PbdiffPcP propagation, leading to tpp — tppaiffPep 
satisfying seismic observations (Fig. 3d). Additionally, the BML seismic 
structure generates waveforms comparable to InSight’s seismic record 
(Supplementary Information Section 4). 

Shear quality factor profiles, Q,(r), can be extracted from our 
inversion output. Similar to viscosity, Q, depends on temperature 
and pressure, and is a function of frequency, w (refs. 23,52,53): 
QC, w) « [n(r)/w]"7, where the power-law exponent a, = 0.1-0.3 
expresses the frequency dependence of Q, (refs. 23,35). The BML Q, 
profiles at 1 Hz (Fig. 4b) have low values in the partially molten mantle. 
Elsewhere, Q, is relatively large for both inversion outputs, indicating 
weak seismic attenuation in the solid mantle. In this region, BML 
models exhibit a steeper increase in Q, (Fig. 4a—b). Deep mantle seismic 
Q,,can also be inferred using the relative attenuation between S and 
ScS phases to infer Q,,s, the apparent Q, along an S wave reflected in 
the lowermost mantle. The approach relies on the amplitude ratios for 
ScS and S phases (Methods and Supplementary Information Section 5), 
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yielding Qes = 1,500 + 500. This contrasts with the relatively attenuat- 
ing Martian mantle at longer tidal periods: Q, = 95 + 10. While BML and 
non-BML models satisfy the constraints on Q, (Fig. 4c), several non-BML 
models with large a, fall outside the range Q, = 95 + 10 while matching 
Qses = 1,500 + 500 (Fig. 4c). This reflects the difficulty of satisfying 
geodetic and seismic constraints across a broad range of time scales. 
BML models are less affected because the partially molten layer accom- 
modates tidal dissipation’. 

An additional observational constraint can be gained by considering 
the gradient of the apparent Q,,in the solid convecting mantle, RG, with 
Ro< indicating a depth decreasing Q,, and vice versa. We estimated 
Ro from the travel times and Qs«s for events sampling different mantle 
depths (Methods), leading to Ro=1.9 + 0.6. Non-BML Ro values are far 
below this acceptable range, consistent with the weak sensitivity of 
Q, in these models (Fig. 4a). In contrast, BML models have higher Rọ 
values, with a large fraction matching the acceptable range, consistent 
with the increase of Q, in the mantle (Fig. 4b). Most models that satisfy 
both Rý and Q,,; constraints have a, = 0.2 — 0.25. These results clearly 
favour BML models, and provide information on the so far poorly con- 
strained frequency dependence of Q,,in Mars’s solid mantle. 

The revised smaller core for BML models is around 5-8% denser than 
previous seismic estimates, implying a reduced light element content 
compared to refs. 1,17 (Methods). Indeed, Mars’s core can be explained 
within cosmochemical bounds by a combination of 17 wt% of sulfur and 
2.9 wt% of oxygen in addition to iron, or bya more complex mixture of 
sulfur, oxygen, hydrogen and carbon in various smaller proportions 
(Methods and Supplementary Information Sections 7 and 8). 

Mars’s magnetic activity went extinct 3.8-4.0 Gyr ago™*°°. Mag- 
netic field production via a thermally driven dynamo action requires 
efficient convective motion in the metallic core, implying core heat 
loss at the CMB. However, the BML heating and heat buffering effect 
prevents core cooling (Fig. 1h). While this suggests that a thermally 
driven dynamo cannot operate on Mars, alternative external sources 
may power an early dynamo. These include early core super-heating 
following core-mantle segregation”, CMB heat flux enhancement due 
to the overturn that led to the formation of the BML, dynamo episodes 
caused by late giant impacts*®, or core elliptical instabilities created by 
satellites orbiting Mars in retrograde fashion”. These mechanisms, 
alone or combined, may have allowed for a Martian dynamo lasting 
for several hundred million years’. 

The presence of a BML on Mars leading to the coexistence of fully 
and partially molten layers above the metallic core requires a reanaly- 
sis of geophysical data, in particular InSight’s seismic record for deep 
phases arrival windows, which may reveal a zoology of seismic phases 
interacting with lowermost mantle discontinuities. 
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Methods 


Inversion results 
Extended Data Table 1 summarizes the inversion results for BML and 
non-BML sets. 

The main and supplementary sets of BML inversions (see Supple- 
mentary Information Section 3 for further details) show compara- 
ble statistics (see columns 4 and 5 in Extended Data Table 1), with the 
exception of the P-wave velocity jump across the CMB of opposite sign 
between the two sets. The velocity jump is negative for the main set and 
essentially positive for the supplementary set, as a result of the smaller 
V,in the liquid silicate layer due to the different equation of state (EoS) 
for the liquid silicates used in this set. Additionally, the models from 
the supplementary set are in good agreement with the seismic data 
and similar to the main BML inversion set (as shown in Fig. 3). 


Attenuation at different frequencies from inversion outputs 
Our two sets of inversion outputs can be used to determine the 
lowermost-to-uppermost viscosity ratio in the solid convecting 
mantle (that is, the mantle regions excluding the lithosphere and the 
BML when present). As mentioned previously, the differences in V* 
imply that the convecting mantle becomes considerably more viscous 
for BML models because pressure and temperature gradients in this 
region are similar for the output models of both inversion sets, but V* 
is larger for BML models. Indeed, neglecting the temperature effect 
(due to the modest adiabatic increase of approximately 130-150 K) 
compared to that of pressure in the convecting mantle, and assuming 
a pressure increase in AP of around 15 GPa throughout the convecting 
mantle, equation (1) implies a lowermost-to-uppermost solid mantle 
viscosity ratio, R, ~ exp[APV*/(RT,,)], leading to R, = 50 for non-BML 
models and R, = 5,000 for BML models. 

The distinct rheology between the homogeneous mantles and BML 
mantles may not be straightforward to explain because of the remaining 
trade-offs between poorly constrained quantities that can affect the 
rheological behaviour (that is, water content, grain size, major element 
composition). For example, an intrinsically more sluggish mantle can 
be explained by a drier mantle*” but also by larger grain sizes, which 
are not explicitly modelled in our inversions. 

As for viscosity, our inversion output can be used to deduce shear 
quality factor, Q, across the planet's radius, r. Q, depends on tempera- 
ture and pressure, and is a function of frequency” >’; 

1 (E*+P(r)V*) "4 
Q, (0,00) = areo ari ) : Q) 
where the power-law exponent a, modulates the dependence of 
attenuation on frequency (w) and Q, is a constant adjusted to match 
constraints on Mars’s global degree-two shear attenuation, Q,, at the 
frequency of the main Phobos tide (5 h 55 min)”. We adjusted Qù, to 
obtain Q, inthe range Q, = 95 + 10 (refs. 3,22,34,35). Qo implicitly con- 
tains information on several rheological parameters (that is, reference 
grain size and exponent or relaxation time scale). While more explicit 
models could be considered for Q,(r, w) (for example, Andrade or 
Burgers“) they would display dependencies on P, T and øo (ref. 35) 
as equation (2) does. 

Because the value of the power-law exponent remains debated, 
with estimates ranging between approximately 0.1and approximately 
0.3 (see refs. 23,35 and references therein), four values of a, were 
considered. 

Inthe solid mantle below the lithosphere, BML models exhibit steeper 
increase in Q, values, with ratios of lowermost-to-uppermost (solid) 
mantle attenuation of Ry = 2-12, while non-BML models are associated 
with much smaller values: Rg = 1-3 (Fig. 4a,b). This difference results 
from the larger values of effective activation volume associated with 
the BML inversion outputs. Using equations (1) and (2) one can relate 


the ratios of lowermost-to-uppermost solid mantle viscosities (R,,) and 
shear quality factors: Rg = Rn! Hence, using the previously estimated 
values of 50 (non-BML) and 5,000 (BML) for R,, and for the range of a, 
considered, this translates to Ry = 2-3 for non-BML and Ry = 4-13 for 
BML inversion outputs. These estimates are close to those displayed 
in Fig. 4a,b. These results can be compared to observational constraints 
on seismic attenuation within Mars’s lithosphere and mantle inferred 
from InSight’s seismic waveforms”. 

To gain information on seismic Q in the deep mantle, one can con- 
sider the relative attenuation between S and ScS phases for suitable 
events to infer Q,.<, the apparent Q value along an S wave reflected in 
the lowermost Martian mantle, without being sensitive to the source 
cut-off. The latter is proposed to be relatively low for quakes originat- 
ing from the Cerberus Fossae region®, which may therefore affect any 
attenuation measurements. This determination relies on the relative 
amplitudes A,., and A, for ScS and S phases, respectively, recorded by 
the seismometer (see following section for details). We conducted 
this analysis using the S1222a event that has the best signal-to-noise 
ratio and an estimated epicentral distance of 37°, for which S and ScS 
phases were detected, yielding Q.,, = 1,500 + 500 associated with this 
event (see section below). 


Qestimation from direct arrivals 
Following ref. 66 we write the amplitude spectrum (A,) of a seismic 
phase ø (for example, S or ScS) as: 


Ag (@) = S(@) Spexpl-(Wty)/(2Q,)1, (3) 
where ©, is a constant and represents the radiation pattern and geo- 


metric spreading of the g phase. Then, using the above expression we 
compute the amplitude spectrum ratio between ScS and S as: 


_ Ages(@) _ Dees _O| tss _ ts 
Rscs-s As(@) Ø; a») 7 Ge al (4) 


This equation expresses the amplitude ratio in theory. Note that the 
ScS is due to an S-wave reflection occurring at the base of the partially 
molten BML, below which the S-wave velocity is zero. For an incident 
Sh wave (like in SO185a), this reflection is a total reflection where all 
the incident Sh-wave energy is reflected as an Sh wave, and the equa- 
tion above holds as is. However, for an incident Sv wave (like in S1222a, 
see Supplementary Information Section 1.2), only a fraction of the 
incident energy is transmitted to the reflected Sv-wave energy, while 
the complementary incident energy will be transferred to the reflected 
Pwave orto the P-wave energy in the fully molten layer below the reflec- 
tor. To characterize the fraction of the incident Sv-wave energy being 
reflected as Sv waves, we compute the amplitude ratio between the 
reflected and the incident Sv waves (that is, using the analytical expres- 
sion for the reflection coefficient, equation 5.40 in chapter 5 in ref. 67) 
at the base of the partially molten BML using typical values of param- 
eters listed in Extended Data Table 2. In addition to the values listed 
in Extended Data Table 2, the base of the BML is considered as a sharp 
boundary and V, below the boundary is O m s”. To avoid numerical 
instabilities in using equation 5.40 of ref. 67 we replaced the zero value 
by10°m s “for V, below the BML. We tested other values (for example, 
10 *mstand10°ms°), which resulted in similar reflection coefficient 
values compared to those obtained with V, = 10° ms“. We found that 
the reflection coefficient is always larger than 98% for the incident 
angle ranging from 0 degrees (vertical incidence) to 10 degrees. Thus, 
we can reasonably neglect these coefficients in equation (4) above. 

The estimation of Q,,, using the equation above requires the know- 
ledge of the travel time of ScS and S phases, Q,, and Ag,</As. For $1222a 
and S0185a, we picked the S based on the Marsquake Service catalogue 
and picked the ScS using waveform matching and polarization analysis’. 
We conducted our picking independently, and our S0185a ScS pick 


agrees with those in ref. 14. For each event, we used a 20-second-long 
time window for the two phases after the corresponding arrival times 
and compute the spectrum ratio in the frequency domain (A,,</As). 
For SO185a, we used the transverse component, as the ScS candidate 
is obvious on the transverse component; we adopted the S1222a radial 
component data because the ScS candidate is more pronounced onthe 
radial component. We used Q; = 4,000 from the direct-arrival Q estima- 
tion (Supplementary Information Section 6). Analysing this informa- 
tion and the equation above, our estimates of Q,,. range between 1,000 
and 2,000. Note that this estimation is essentially model-independent: 
even though the absolute arrival time t, is model-dependent we checked 
that changing its value by 100 seconds did not affect the lower and 
upper bounds of 1,000 and 2,000 for Q,,s. All other parameters that 
enter equation (4) are based on seismic observations. 

Most of the BML and non-BML models satisfy the constraint on Q}. 
This is not surprising given that Q, is adjusted to match the required 
range for Q, (Fig. 4c). However, several non-BML models with relatively 
large a, values do not fall within the acceptable range for Q, even if 
they agree with the Q,., range (Fig. 4c). BML models are less affected 
because the presence of a deep partially molten layer concentrates 
the dissipation in this region at Phobos tidal frequency’. Note that a 
partially molten and strongly attenuating layer with a typical thickness 
of 100 km or less would correspond to only about 10 wavelengths or 
less for a 2-second-period body-wave, which would therefore not be 
appreciably damped by the partially molten structure. Furthermore, 
the presence of the soft basal layer favours a smaller core radius to 
satisfy constraints on the real part of the tidal Love number, k,, as 
demonstrated for the Moon®. 

Since a large fraction of BML and non-BML models fall within the 
acceptable Q, and Q,,, ranges, it remains difficult to discriminate 
between the two sets of models because Q,,; values remain dominated 
by large (shallow) values of Q, along the ScS paths. 

Nevertheless, an additional observational constraint can be gained 
by considering the gradient of the apparent shear attenuation in the 
solid (that is, convecting) mantle. As seen above, such gradient (Ro) is 
considerably larger for BML models (Fig. 4a,b). If one considers another 
event witha larger epicentral distance to compute Q,,<, the correspond- 
ing S and ScS raypaths will sample deeper regions of the mantle 
(Fig. 1f,1), which may result in different Q,., values compared to those 
associated with S1222a, especially for BML output models that exhibit 
a larger increase of Q,, with depth (thatis, larger Rg values). We applied 
this approach to event S0185a that produced observable ScS and S 
phases and has a larger epicentral distance than S1222a (55 degrees for 
S0185a versus 37 degrees for S1222a), requiring an S-raypath diving 
into deeper regions of the mantle. Using these two events, the ampli- 
tude ratio Rý between (Ages/As) 777" and (Ages/As)>"8*" provides a 
constraint on the gradient of the shear quality factor in Mars’s solid 
mantle. A ratio Rý smaller than unity indicates a decreasing Q, with 
depth, and vice versa. Due to the integrated nature of Q, and Q,,, there 
is no direct correspondence between Rý and Ro. However, a relative 
correspondence exists: an increase in Rg yields an increase in RG, with 
the converse also being true. One can estimate Rgas a function of travel 
time and apparent Q for ScS and S phases for events $1222a and SO185a. 
Making the reasonable assumption, as justified above, thatØses /Øs <1, 
and using equation (4), yields the following expression for relative 
ratios between the two events: 


Ro= exp[@(Atg- Atg,s)/2], (5) 


where At} (with g being either the phase S or ScS) is the difference 
of the ratios of the travel time to the apparent attenuation, 
t,/Q,, between events S0185a and S1222a, which are computed along 
each corresponding raypath. Note that here the computation of the 
relative ratios is purely based on the models. The observational range 
for Ro was determined by fitting independently R,,s., values obtained 


in the frequency range 0.2-1.0 Hz using linear expressions. To deter- 
mine the amplitude ratio Rs«s-s we used the same S and ScS time win- 
dows fromthe ScS waveform matching of S0185a and S1222a. We taper 
the time windows using a cosine function and then convert the tapered 
waveforms to the frequency domain to obtain the amplitude spectra. 
For each event, we compute the amplitude ratio between the ScS 
and S using the corresponding spectra. A few anomalously high peaks, 
with amplitudes approximately 10 times larger than the Rss- values 
measured for the rest of the population (two peaks for S1222a and one 
for SO185a), were found at certain frequencies. These peaks are most 
likely due to noise contamination and were therefore removed before 
performing the linear fitting. Then, we computed the L, error norm on 
the fit for both S1222a and S0185a, leading to Rọ =1.9+0.6. 


Constraints on core composition 

The higher core density inferred for BML models implies a reduced 
content in light elements compared to what was recently suggested in 
refs. 1,17. The most likely light elements in the core of Mars are, inorder 
of abundance, sulfur (S), oxygen (O), carbon (C) and hydrogen (H)””. 
They are soluble in iron during core formation and abundant enough 
to affect the elastic properties of the core. To deduce the composition 
of the core, we first defined the subset of possible Fe-O-S-C-H alloy 
compositions that are consistent with the nickel and cobalt 
trace-element composition of the Martian mantle” by massively sam- 
pling core formation models” (Supplementary Information Section 7). 
Then, we searched for those compositions whose EoS best agrees with 
the isentropic EoS parameters of the Martian core inferred in this study 
(density (Ocyp), isentropic bulk modulus (K;) and its pressure derivative 
(K’) at CMB conditions, see Supplementary Information Section 8). 
The EoS for the liquid alloy was built using the most recent experimen- 
tal data relevant for the composition, pressure and temperature con- 
ditions of the Martian core. This is the most suitable approach to 
account for conditions prevailing in Mars’s core that cannot yet be 
captured by alternative ab initio approaches (Supplementary Informa- 
tion Section 8). We focused this analysis to BML models only as non-BML 
models are discussed in previous recent works”. To deduce the com- 
position of the core we only retain models for which the EoS parameters 
agree with those of the alloy considered here within 0.5%, 1% and 15% 
respectively (Supplementary Information Section 8). For all models 
the average core density can be matched by a Fe-0-S-C-H alloy using 
different proportion of O,S,C and H, but only a small fraction of those 
models agree with the acoustic velocity of the alloy (see Supplementary 
Information Section 8). The revised core size we inferred (R, = 1,650 + 
20 km) lies at the lower end of the InSight pre-mission range”, is about 
170 km smaller than recent seismically derived estimates and yields 
an average core density p. = 6.5 g cm”, which is about 5-8% denser 
than previous seismic estimates and can be explained by an iron alloy 
with a smaller amount of light elements. Contrary to previous recent 
studies”, the core of Mars can be explained by metal-silicate data 
within cosmochemical plausible bounds (for sulfur) by a mixture of Fe 
and17 wt% of Sand 2.9 wt% of O, or bya more complex mixture of S, O, 
H and C in various smaller proportions. Fe-S-O-C assemblages are 
most likely composed of 14 wt% of S, 3 wt% of O and 1 wt% of C. The 
addition of H (inamounts of 0.15 wt%, or even 0.5 wt%, that is, beyond 
the maximum bound allowed by experimental petrology and cosmo- 
chemistry) further decreases the fractions of other light elements. 
However, the effect of adding H to Fe-S-O-Cis, on average, less impor- 
tant (that is, causing a 1 wt% decrease or less of the fractions of S, O and 
C) because the small quantities of H permitted here do not appreciably 
affect the proportions of other light elements (Supplementary Infor- 
mation Section 8). 

Extended Data Table 3 summarizes the statistics of the results 
(also displayed in Supplementary Fig. 11), with the addition of a case not 
shown in in Supplementary Fig. 11, which corresponds to a Fe-O-S-H 
assemblage where the maximum upper bound for H was set to 0.15 wt% 
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(as determined in Supplementary Information Section 7.3) instead of 
0.5 wt% in the fourth column. Fe-O-S-C assemblages are composed of 
15+2wt% of S,3 +1wt% of O and1+1wt% of C. For the more complex 
assemblage containing 0.3 + 0.1 wt% of Hin addition to the three other 
light elements, their proportions decrease slightly by 1-2%. Interest- 
ingly, an Fe-O-S-H that contains less than 0.15 wt% of H results, on 
average, in similar composition for the other light elements: 14 + 2 wt% 
of S,3 +2 wt% of O and 1 +1 wt% of C. These fractions are also similar 
to those of the simpler assemblage Fe-O-S-C. This is because the 
small amounts of H allowed do not significantly impact the other light 
element contents. 
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Extended Data Table 1| Summary of the inversion results with a homogeneous mantle and with a BML (best 1,000 output 
models for each inversion set) 


Parameter Meaning Non-BML BML (main set) BML (Suppl. set) Units 
E* Mantle effective activation energy 320+70 120+15 120+10 kJ/mol 
v* Mantle effective activation volume 4.642 7.32 8.141 cm? /mol 
No Mantle reference viscosity 1072+0.4 1070:2= 0:4 1040:9205 Pa s 
ka BML thermal conductivity - 51 62 W m7? K! 
Ling Initial uppermost mantle temperature 1750+30 1820+60 181070 K 
Teo Initial CMB temperature 2230+80 2240+75 2240+100 K 
Tp Potential temperature 1790+40 151020 1530+20 K 
Tm Uppermost mantle temperature 1880+50 1540+20 1610+130 K 
Ta Average BML temperature - 2390+90 2310+150 K 
Te CMB temperature 2090+40 2760+150 2570+150 K 
Der Crustal thickness 71+1 6742 6743 km 
Diu Lithospheric and upper TBL thickness 600+70 260+20 250+10 km 
Da BML thickness - 165+20 140+20 km 
Re Core radius 1820+10 165020 1660+20 km 
Pe Pressure at CMB 19+0.2 22+0.3 21.5+0.3 GPa 
Vp. eis Vp in the mantle at CMB 9.37+0.09 5.59+0.04 4.86£0.03 km/s 
Bic Vp in the core at CMB 4.87+0.12 4.99+0.10 5.06+0.09 km/s 
AVpoun Vp —Vp at CMB -4.50+0.16 -0.59+0.11 0.200.080 km/s 
PCMB Core density at CMB 5870+40 6240+60 6240+60 kg/m? 
Ks Core isentropic bulk modulus at CMB 1407 16010 1606 GPa 
K5 Pressure derivative of Ks 5.1+0.8 6.2+0.5 5.7+0.7 = 
Pc Average core density 6130+30 6470+60 6470+60 kg/m? 
Quantities listed correspond to constant or to present-day values unless specified otherwise. The subscripts ‘0’ indicate initial values for quantities that vary in time. 
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Extended Data Table 2 | Typical P- and S-waves velocity and density values used to estimate the S-wave reflection coefficient 
for models with a BML 


Layer Vp [km/s] Vg [km/s] Density [kg/m3] 
Above the BML base 7.1 1.4 4800 
Below the BML base 5.5 0 4800 


Typical values of V,, V, and density above and below the BML used to calculate the amplitude ratio between the reflected and the incident Sv-waves at the base of the partially molten BML. 


Extended Data Table 3 | Statistics for the core compositions that match seismic data for four assemblages of light elements 
considered (average value and 1-sigma range) 


Element [wt.%] | Fe-O-S Fe-O-S-C  Fe-O-S-H  Fe-O-S-H (H<0.15 wt%) 
S 17.2+0.9 14.6+1.7 13.041.8 14.1+1.8 
O 2.9+0.5 2.7+0.5 2.3+0.5 2.60.5 
C - 1.10.5 1.30.6 1.30.5 
H - - 0.30.1 0.10.0 


The last column corresponds to a Fe-O-S-H with a H content below the absolute maximum bound of 0.15 wt%. 


